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RENORMALIZATION  OF  AN  INVERSE  SCATTERING  THEORY 
FOR  INHOMOGENEOUS  DIELECTRICS 


INTRODUCTION 

Inverse  scattering  theory  seeks  to  reconstruct  the  unknown  physical  properties  of  an  object  from 
information  contained  in  the  scattering  data.  In  this  paper  we  consider  the  reconstruction  of  the  dielec¬ 
tric  profile  of  an  inhomogeneous  slab  from  the  electromagnetic  reflection  coefficient.  Approximate 
solutions  to  this  problem  are  often  obtained  by  considering  a  limited  range  of  frequencies  or 
wavelengths.  For  example,  the  wave  equation  for  the  electric  field  can  be  converted  into  an  integral 
equation  which  is  then  solved  by  iteration.  In  the  first  iteration  the  Born  approximation  is  the  solution 
obtained  by  assuming  small  phase  shifts  in  the  scattered  field.  The  radius  of  convergence  for  this 
approximation  is  limited  to  very  small  values  of  the  wave-number,  so  that  inversion  methods  based  on 
this  approximation  will  have  a  limited  range  of  applicability. 

In  this  paper  we  investigate  a  method  to  increase  the  radius  of  convergence  of  approximate  solu¬ 
tions  of  inverse  scattering  problems  by  using  renormalization.  We  base  our  inversion  method  on  the 
exact  theory  that  has  been  developed  by  Kay  and  Moses  [1,2,3]  from  the  mathematical  investigations  of 
Gelfand,  Levitan  and  Marchenko  [4,5].  The  search  for  approximate  solutions  to  be  used  for  the  elec¬ 
tric  field  within  inhomogeneous  dielectric  layers  led  us  to  consider  the  method  of  multiple  scales  [6], 

By  considering  the  continuity  of  the  electric  field  and  its  derivative  at  the  boundary  of  the  slab 
and  by  using  the  renormalized  electric  field,  a  rapidly  convergent  solution  is  obtained  that  is  equivalent 
to  a  second-order  perturbation  treatment  of  the  exact  inverse  scattering  theory.  In  the  lowest  order  the 
renormalized  electric  field  is  equivalent  to  the  WKB  approximation  [7].  This  theory  is  demonstrated  by 
numerical  examples  using  data  obtained  from  a  numerical  simulation  of  scattering  from  Gaussian  and 
parabolic  dielectric  profiles  and  homogeneous  slabs. 

INVERSION  THEORY 


We  consider  the  idealized  physical  model  shown  in  Fig.  1.  A  time-harmonic  plane-polarized  elec¬ 
tromagnetic  wave  of  wave  number  k  —  2 n/\ 


Em  ( x,k )  -  eikx  ,  i  s  V=T,  x  <  -1/2  , 


(1) 


is  normally  incident  from  the  left  upon  an  inhomogeneous  dielectric  slab  of  thickness  L .  The  permit¬ 
tivity  relative  to  free  space,  er,  is  a  function  of  the  geometric  distance  x;  we  first  demonstrate  our 
method  with  e,  defined  by  the  truncated  Gaussian  function 


«,  (x) 


Ae-*1 ,  |x|  <  1/2 

,1.0  ,  I x i  >  L/2  ’ 


(2) 


where  b  is  chosen  to  ensure  continuity  of  «,  at  |x|  —  1/2.  The  wave  amplitude  E(x,k)  satisfies  the 
scalar  Helmholtz  equation 
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fE(*'k)  +  k2  e,(x )  E(x,k)  -  0. 


(3) 


The  reflected  wave  amplitude  is 

Esc(x,k)  -  r(k)  e~ikx  ,  x  <  -LI 2  (4) 

where  r(k)  is  the  complex  reflection  coefficient. 

The  electromagnetic  path  length  s  in  the  slab  is 

s  =  s(x,Ar)  -  f_L/2  k  y/er(x')  dx\  (5) 

which  is  the  Liouville  transformation  between  geometric  x-space  and  path-length  s -space  [8].  The 
inverse  problem  that  we  consider  is  based  upon  the  Helmholtz  Eq.  (3)  with  an  arbitrary  relative  permit¬ 
tivity  er(x).  If  «,(x)  is  a  monotonic  function  of  x  and  independent  of  k,  differential  Eq.  (3)  can  be 
transformed  into  a  Schroedinger-type  equation  in  s-space  whose  potential  function  q(s)  can  be  deter¬ 
mined  uniquely  if  r(k)  is  an  analytic  function  of  k  [1,2].  This  theory  has  been  well-developed, 
although  the  transformation  to  obtain  a  unique  «,(x)  from  q(s)  depends  on  the  functional  form  of  er, 
as  is  evident  from  Eq.  (5),  and  is,  in  general,  nontrivial. 

The  inversion  method  will,  perforce,  first  reconstruct  «,  in  s-space  and  will  then  make  a  transfor¬ 
mation  to  x-space,  based  on  a  priori  knowledge  of  e,(x).  For  notational  convenience,  we  indicate  the 
different  functional  forms  of  er  in  these  two  spaces  by  reference  to  the  argument  x  or  s.  The  Lorentz 
model  for  the  permittivity  of  a  collection  of  free  electrons  is  a  simpler  special  case,  where 

er(x,k)  -  1  -  q(x)/k2  (6) 

and  where  the  profile  function  q(x),  which  is  the  electromagnetic  analogue  of  the  quantum-mechanical 
potential  function,  is  proportional  to  the  electron  density.  Substitution  of  Eq.  (6)  into  the  Helmholtz 
Eq.  (3)  leads  directly  to  a  Schroedinger-type  equation  in  x-space. 

Referring  to  the  Helmholtz  Eq.  (3),  if  we  could  assume  that 

i/»(s,k)  -  |  €r(x(s))  J  1/4  E(s,k)  (7) 

then  it  would  be  possible  to  transform  Eq.  (3)  into  a  Schroedinger  equation  of  the  form 

--  +  [  k2  -  q(s )  ]  1 l/(s,k)  -  0  (8) 

and  where  the  function  q(s)  is  to  be  determined.  This  assumes  we  know  the  functional  form  of  «,  and 
can  evaluate  the  integral  in  Eq.  (5)  to  obtain  x(s)  explicitly;  for  an  inverse  problem,  this  a  priori  infor¬ 
mation  would  not,  in  general,  be  available. 

The  Fourier  transform  V(s,t)  of  >J/(s,k)  satisfies  the  time-dependent  wave  equation 

-|4'F(s,f)--^-'I'(s,/)-(7(s)'lr(s,t)-0  (9) 

os*  dr 

where  t  is  the  time  variable  with  the  velocity  of  light  c  =  1.  In  free-space,  x  <  -L/2,  the  incident 
plane  wave,  Eq.  (1),  and  the  incident  wave  amplitude  for  Eq.  (8),  are  identical;  the  incident 

plane-wave  for  the  time-dependent  Eq.  (9)  is  represented  by  a  unit  impulse. 


Vin(s,t)-8(s-  f) 


(10) 


which  produces  the  reflected  transient,  or  characteristic  function, 

R(s  +  t)~  4~  f~  '•<*>  e~ikis  +  ,)dk. 
2 it  J~°° 


(11) 


Due  to  causality,  we  must  have 

R(z)  -  0 , 2  <  -L/2  (12) 

i.e.,  a  reflected  transient  is  not  produced  until  the  incident  pulse  has  interacted  with  the  inhomogeneous 
medium.  It  is  possible  to  relate  the  wave  amplitude  Wisj)  in  the  inhomogeneous  slab  with  the  wave 
amplitude  ¥(,($,/)  in  the  free-space  region  by  the  linear  transformation  [1] 

¥($,/)  -  Vo(s,t)  4-  K(s,z)  'Po (z,t)  dz ,  (13) 

where 

-8 (s  -  t)  +  R(s  +  t).  (14) 

From  physical  considerations  we  know  the  V(sj)  is  a  right-moving  transient,  so  that 

'V(sj)  -  0  ,  s  >  t.  (15) 

Thus  ATs,r)  =•  0  for  t  >  s.  Substituting  the  expression  (14)  in  (13)  and  using  Eq.  (12)  and  (15) 
yields  the  integral  equation 

K(sj)  4 -R(s  +  t)  +  /*  AT(s,z)  R (z  +  t)  dz  -  0,  (16) 

which  is  Kay’s  version  of  the  Gelfand-Levitan-Marchenko  integral  equation  of  inverse  scattering  [31. 
Substitution  of  the  representation  (13)  in  the  wave  Eq.  (9)  shows  that  the  function  K(s,t)  satisfies  the 
same  differential  equation  as  the  wave  amplitude  ¥(s,t)  if  the  following  conditions  are  imposed, 

K(.s,-s)  —  0,  and  2~  K(s,s)  —  ?(s).  (17) 

as 

If  the  integral  Eq.  (16)  can  be  solved  for  the  function  K(s,t)  then  the  second  condition  in  Eq.  (17) 
yields  the  exact  solution  to  this  inverse  scattering  problem  [9].  Moses  [2]  has  obtained  the  second- 
order  regular  perturbation  solution  for  q(s): 

q(s)  "  “2*  *(2s)  +4  [  *(2s)  1  (18) 

This  approximate  solution  of  the  inverse  problem  for  a  Schroedinger  equation  can  provide  a  starting 
point  to  reconstruct  «,(s)  [10,11]. 

Here  we  go  directly  to  the  inverse  problem  by  considering  the  electric  fields  within  the  slab.  We 
note  that  by  using  the  Liouville  transformation  (5)  the  differential  Eq.  (3)  for  the  electric  field  within 
the  slab  can  be  expressed  in  terms  of  the  electromagnetic  path  length  s  as 

-  E[s)-  4-  g(x)  4-  E(s)  -  0,  (19) 

ds 2  ds 


3 


g(x )  =  -J£  er(x ) 


-3/2  d 6,0t) 
dx 


If  the  function  g(x)  is  small  in  magnitude, 

max  |g(x)|  =  y  <  1,  Ul  <  1/2  (21) 

then  we  can  apply  a  regular  perturbation  expansion  in  powers  of  y  to  eq.  (19).  However,  such  an 
expansion  would  not  be  uniformly  valid  due  to  the  presence  of  secular  terms  in  the  higher-order 
approximations.  Moreover,  the  explicit  functional  form  of  g(x)  is  unknown  in  the  inverse  problem. 

In  performing  the  transformation  of  the  Helmholtz  Eq.  (3)  into  the  Schroedinger  Eq.  (8),  it  was 
implicitly  assumed  that  x  can  be  determined  from  eq.  (5).  The  method  of  multiple  scales  is  used  to 
provide  an  effective  technique  for  summing  perturbation  sequences  when  the  explicit  functional  form 
of  €,(x)  is  not  known.  Using  this  method,  we  seek  an  approximate  solution  for  the  electric  field  within 
the  slab  in  the  form 

E(s(x,k),x)  «■  E0(s,x)  +  yE\(s,x)  +  y2£2(s,x)  +  . . .  ,  (22) 


where  5  and  x  are  regarded  as  the  fast  and  slow  variables  respectively.  The  differential  equation  for 
E(s,x)  can  now  be  written  through  the  first  order  in  the  perturbation  parameter  y. 

yO  .  £0(s  x)  +  £0(s,x)  =»  0  (23) 


ds 2  u 
d2£!(s,x) 
ds2 


+  £[(s,x)  -  -£(;t) 


dE0(s,x ) 


ft ,(x)  dsdx 


E0(s,x). 


The  general  solution  of  Eq.  (23)  is  given  by 

E0(s(x,k),x)  -  /40(x)  e*  +  B0(x )  e~i:.  (24) 

Since  we  are  only  interested  in  right-traveling  waves,  we  set  Bo(x)  —  0.  The  zero-order  solution  is  cou¬ 
pled  with  the  first-order  equation  for  the  electric  field  amplitude.  This  coupling  leads  to  secular  growth 
in  f^s,.*),  which  violates  energy  conservation.  In  order  to  eliminate  this  growth  the  functional  form 
of  ,40(*)  must  be  determined  from 


y [  fr(x) 


-3/2  der 

~dx 


.  ,  „  f  .  ,  l-i/2  BAq  8£0 

+  2  [  | 


which  has  the  solution 


A0(x ) 


[  er(x) 


so  that  the  renormalized  form  for  E0(s,x)  becomes 

E0(s,x)  -  e"[  fr(x)  1  _1/4,  s  -  s(x,k).  (27) 

The  functional  form  (27)  was  used  to  transform  the  original  Helmholtz  Eq.  (3)  into  a  Schroedinger 
equation.  The  multiple  scale  technique  demonstrates  that  this  functional  form  arises  from  denying  sec¬ 
ular  growth  in  the  first-order  correction  for  the  electric  field  amplitude.  This  process  of  summing  per¬ 
turbation  expansions  to  make  them  more  uniformly  valid  is  called  renormalization  [12]. 


Equation  (19)  can  be  converted  into  an  integral  equation  by  using  an  appropriate  Green’s  func¬ 
tion.  Consideration  of  the  asymptotic  form  of  the  wave  solution  gives  an  exact  expression  for  the 
reflection  coefficient, 

r(k)  -  f  [€,(x)  l-3/2  -j-  (  £(s(x,fc),x)  I  ds.  (28) 

2  J-Ll2  dx  ds  '  t 

Substituting  the  renormalized  solution  (27)  into  Eq.  (28)  and  application  of  the  Fourier  transform  gives 
the  approximate  solution  for  «,(s)  in  s- space 

e,(s)  -  expj  -4 J  ^  R  (z)  dz J  (29) 

where  the  R(,z)  term  in  the  exponent  represents  the  Fourier  transform  of  the  scattering  data,  as  dis¬ 
cussed  in  the  next  section. 

A  complete  inversion  method  should  also  be  able  to  obtain  «,  as  a  function  of  the  geometric  dis¬ 
tance  x,  as  discussed  by  Kay  [1].  The  Born  approximation  relates  «,  in  x-space  to  the  Fourier 
transform  of  r(k).  It  is  possible  with  our  inversion  method  to  transform  from  s-space  to  x-space  if  the 
slab  thickness  L  is  known.  The  WKB  approximation  corresponds  to  the  first  term  of  the  perturbation 
solution  in  k -space. 


Since  we  want  to  relate  a  discrete  set  of  data  points  for  e,(s)  to  the  corresponding  points  within 
the  slab  in  x-space,  we  consider  the  inhomogeneous  region  to  composed  of  thin  dielectric  s:  ■'  ,  of 

thickness  A,  so  that  x,  =  -y-  +  (/  —  1)  A,  /  —  1 . N.  The  number  of  data  points  N  can  be  deter¬ 

mined  since  we  also  know  the  beginning  and  end  points  of  the  slab  in  both  x-  and  5 -spaces;  namely. 

r  L/2  _ ^ 

x(  -  -LI 2  corresponds  to  s  -  0  and  xN  -  L/2  corresponds  to  s  -  k  J_L/2  y/*,(x')  dx\  as  shown  in 

Fig.  3.  In  this  example,  N  can  be  chosen  conveniently  to  divide  the  region  in  s-space  where  «,  differs 
from  the  permittivity  of  free  space;  in  this  example,  N  is  approximately  50.  In  order  to  reconstruct  the 
dielectric  profile  in  x-space,  we  expand  the  profile  in  a  Maclaurin  series  about  the  center  of  the  slab 

er(x)  -  € T (0)  +  e’(0)  x  +  «,(())  x2/2  +  ...  .  (30) 


The  data  points  for  er(s)  obtained  from  Eq.  (29)  can  be  transformed  to  the  equivalent  data  points  for 
€ r (x)  by  fitting  to  a  least  squares  polynomial  whose  coefficients  can  be  related  to  Eq.  (30). 


GENERATION  OF  SCATTERING  DATA 


The  experimental  reflection  data  are  simulated  by  numerical  integration  of  the  exact  Riccati  equa¬ 
tion  for  the  complex  reflection  coefficient,  r(k).  The  direct  scattering  problem  is  modeled  by  a  single 
slab  whose  permittivity  profile  is  a  Gaussian  function  centered  at  the  midpoint  of  the  inhomogeneous 
region,  as  shown  in  Fig.  1.  The  magnitude  of  the  relative  permittivity  ef  is  matched  to  the  relative  per¬ 
mittivity  of  free  space  at  the  boundaries  of  the  slab,  ix|  -  L/2.  The  reflection  coefficient  at  the  left 
face  of  the  slab  r\k)  -  r(k,-L/ 2),  which  is  used  as  a  data  point  for  inversion,  is  obtained  from  the 
solution  to  the  Riccati  equation  [13] 

■—  (x,fc)  -  — —  |  ef(x)  [  l  +  r  ]2  +  [  1  -  r  ]2  j  ,  |x|  <  y  (31) 


with  the  boundary  condition  r(k,+L/2)  -  0.  The  reflection  coefficient  rik)  is  calculated  at  discrete 
points  x,  through  the  slab  for  each  value  of  wavenumber  k. 


We  use  a  fifth-order  Runge-Kutta  method  with  500  increments  in  x  to  compute  the  complex 
reflection  coefficients,  r (Ac),  at  discrete  wavenumbers,  kp  —  (1  —  p)  8,  p  —  1, . . .  ,M  (in  our  case 
M  —  256  points  were  used  for  k  —  0  to  25,  so  that  the  increment  8  —  0.098).  The  modulus  of  r(k)  is 
shown  in  Fig.  2  for  a  typical  set  of  parameters. 

For  use  in  the  inversion  method,  this  set  of  reflection  data  (31)  is  then  represented  by  a  complex 
Fourier  series, 

a0  m-z  1 

r(kp)  -  —  +  £  am  cos(<t>m)  -ibm  sin(<£m)  +  —  aM  cos (irkp)  (32) 

L  m  -  t  1 

where  <t>m  —  mkpir!p ,  M  -  256  and  P  is  the  half-period  of  the  reflection  coefficient.  We  note  that  our 
inversion  algorithm  assumes  r  (Ac)  is  periodic;  as  shown  in  Fig.  2,  P  -  25,  approximately.  The  real  part 
of  the  reflection  coefficient  is  represented  by  a  truncated  cosine  series  whose  coefficients  are  deter¬ 
mined  by  requiring  that  the  data  points  collocate  with  the  series  at  the  appropriate  wavenumbers;  the 
imaginary  part  of  the  reflection  data  is  represented  by  a  sine  series,  so  the  r*(k)  —  r(—k),  k  real.  This 
continuous  representation  of  the  complex  reflection  coefficient  enables  interpolation  and  accurate 
numerical  determination  of  the  minima  and  maxima  of  r(k).  Moreover,  the  Fourier  representation 
enables  an  analytic  evaluation  of  the  function  R  (s),  which  is  needed  for  the  solution  of  the  integral  Eq. 
(16). 

DISCUSSION  OF  RESULTS 

We  have  demonstrated  an  inversion  method  for  scattering  data  from  a  smooth  dielectric  profile. 
The  scattering  data  of  Fig.  2  were  used  in  our  inversion  method  to  reconstruct  the  permittivity  profile 
of  Fig.  3.  Table  1  compares  the  reconstructed  results  with  several  known  Gaussian  profiles;  the  param¬ 
eters  A,b  were  obtained  by  a  regression  analysis  on  Gaussian  profiles  in  order  to  test  the  accuracy  of 
this  inversion  method  for  increasing  gradients  in  the  dielectric  inhomogeneity.  This  inhomogeneity  is 
measured  by  the  parameter  y,  which  can  be  shown  from  Eq.  (19)  to  be  equal  to  6/2.  As  expected, 
good  agreement  was  obtained  for  smaller  values  of  b  and  poorer  agreement  as  b  was  increased.  The 
correlation  coefficient  p  tests  the  validity  of  the  regression  model;  p  -  1  indicates  that  the  recon¬ 
structed  profile  agrees  perfectly  with  the  model  profile. 

Table  2  tests  the  accuracy  of  our  method  when  applied  to  three  different  profiles:  Gaussian,  pa¬ 
rabola,  and  homogeneous  slab.  Simulated  numerical  scattering  data  were  generated  by  solving  for  each 
of  these  profiles,  as  characterized  by  equivalent  A  and  b  parameters.  The  inversion  method  was  then 
applied  to  the  three  sets  of  data  and  profiles  of  er  were  obtained  in  s-space  and  in  x-space.  By  assum¬ 
ing  that  the  unknown  profile  for  er(x)  can  be  developed  as  a  series  of  the  form  (30),  the  coefficients 
of  the  least-squares  polynomial  can  be  related  to  the  parameters  A ,  b  of  an  assumed  Gaussian,  para¬ 
bolic,  or  homogeneous  slab  model.  A  regression  analysis  with  the  three  model  profiles  was  used  to 
compare  the  profile  parameters  A,  b\  the  correlation  coefficients  p  are  also  shown  in  Table  2.  We  note 
that  the  b  parameters  for  the  parabola  and  the  homogeneous  slab  do  not  ensure  continuity  of  er  at 

Ul  —  y.  For  example,  the  approximation  for  the  Gaussian  profile  is  written 

€r(x)  =  A  (1  -  bx2).  (33) 

The  element  along  the  diagonal  show  results  due  to  "good  guesses",  i.e.,  the  least-square  analysis  com¬ 
pared  the  reconstructed  profile  with  the  same  model  profile  that  was  used  in  the  direct  calculation  of 
the  reflection  data.  The  off-diagonal  elements  shows  results  comparing  "wrong”  regression  models  with 
different  direct  models.  Note  that  this  matrix  is  not  symmetric. 

From  these  results  we  see  that,  to  second  order,  it  is  not  clear  whether  the  Gaussian  or  parabola  is 
the  better  fit;  the  homogeneous  slab  is  clearly  a  poor  fit.  The  correlation  coefficients  for  the  regression 


analysis  for  the  homogeneous  slab  have  not  been  display  since  the  least-squares  fit  is  a  simple  arith¬ 
metic  average,  so  the  correlation  coefficient  cannot  be  defined.  The  qualitative  features  of  the  slab  can 
be  reconstructured  by  the  inversion  algorithm,  but  the  magnitude  of  er  is  too  large  by  20-30%.  This  is 
expected  since  the  homogeneous  slab  presents  a  sharp  discontinuity  that  violates  condition  (21)  for 
convergence  of  the  perturbation  series.  We  conjecture  that  our  inversion  method  can  be  extended  to 
discontinuities  and  steep  gradients  by  boundary-layer  theory. 

The  results  of  Table  1  show  that  the  renormalized  inversion  theory  has  a  larger  radius  of  conver¬ 
gence  than  the  Born  approximation.  Our  inversion  method  also  includes  phase  information  that  is 
neglected  in  the  Born  approximation  [141, 

(kL)  (er(x)  -  1)  <  <  1  .  (34) 

We  note  that  for  L  -  2.0  and  A  =  1.80,  as  used  in  our  examples,  the  Bom  approximation  is  not  valid 
for  wavenumbers  greater  than  0.60  in  Fig.  2.  Profiles  reconstructed  with  the  Born  approximation  over 
a  range  of  k- space  of  0  to  25  were  in  complete  disagreement  with  the  known  profiles.  Renormalization 
has  extended  the  radius  of  convergence  in  Ar-space  by  two  orders  of  magnitude.  The  results  of  Table  2 
suggest  that  the  regression  analysis  can  be  extended  to  higher  orders  to  consider  many  different  model 
profiles.  Extensions  of  the  inversion  theory  to  include  dielectric  absorption  and  scattering  at  oblique 
incidence  would  require  generalizing  the  solution  of  the  integral  Eq.  (16). 
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Input  parameters 

Reconstructed  parameters 

A 

b 

A 

b 

P 

1.50 

0.41 

1.59 

0.41 

0.99 

1.80 

0.59 

1.94 

0.60 

0.97 

4.00 

1.39 

4.70 

1.25 

0.95 

6.00 

1.79 

7.19 

1.55 

0.95 

8.00 

2.08 

9.54 

1.59 

0.94 

10.00 

2.30 

11.54 

1.87 

0.92 

Table  2  —  Comparison  of  Inversion  Results  for  Parameters  A,  b 
for  Gaussian  and  Parabolic  Profiles  and  Homogeneous  Slabs. 
Input  Gaussian  Profile  Parameters:  A  =  1.80,  A  *  0.59,1  -  2.0 


